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ABSTRACT 

We generate stellar distribution functions (DFs) in triaxial haloes in order to examine 
the reliability of slopes T = AlogAf/Alogr inferred by applying mass estimators of the 
form M oc R e cr 2 (i.e. assuming spherical symmetry, where R e and a are luminous effec- 
tive radius and global velocity dispersion, respectively) to two stellar sub-populations 
independently tracing the same gravitational potential. The DFs take the form f(E), 
are dynamically stable, and are generated within triaxial potentials corresponding 
directly to subhaloes formed in cosmological dark-matter-only simulations of Milky 
Way and galaxy cluster haloes. Additionally, we consider the effect of different tracer 
number density profiles (cuspy and cored) on the inferred slopes of mass profiles. For 
the isotropic DFs considered here, we find that halo triaxiality tends to introduce an 
anti-correlation between R e and a when estimated for a variety of viewing angles. 
The net effect is a negligible contribution to the systematic error associated with the 
slope of the mass profile, which continues to be dominated by a bias toward greater 
overestimation of masses for more-concentrated tracer populations. We demonstrate 
that simple mass estimates for two distinct tracer populations can give reliable (and 
cosmologically meaningful) lower limits for T, irrespective of the degree of triaxiality 
or shape of the tracer number density profile. 
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1 INTRODUCTION 

The Milky Way's dwarf spheroidal (dSph) satellites in- 
clude the most dark-matter-dominated galaxies known, with 
dynamical mass-to-light ratios rangin g from order ~ 10 
to several hundreds in solar units (Mateo 1998). This 
makes dSphs objects of prime interest for studying the dis- 
tribution of dark matter in galaxies. dSphs lack atomic 
hydrogen; therefore methods for measuring dSph masses 
must rely on the kinematics of their pressure-supported 
stellar populations. In the past decade, many techniques 
have been developed with the goal of determining the in- 
ternal mass distributions of dSphs: spherical Jean s mod- 
elling dLokas fc Mamonl l200ll; IStrigari et all l200d. I 200S ; 



Koch et al II2007I : IBattaglia et al.ll200Sl ; IWalker et allboOS 



Wolf et al.ll2010l), phase-space modelling (| Wilkinson et al.l 
20021 : Amorisco fc Evansl 2011), the multip le stellar popu- 
lations method l| Walker fc Penarrubial [201lh . the use of the 
virial theorem for sphe rical and constant flattened systems 
l|Agnello fc EvansH2012T ) as well as axisymmetric Jeans mod- 



elling (IHavashi fc Chiball2012f) and Schwarzschild model ing 
i|jardel fc Gebharddl2012l ; iBreddels et al.ll2012T l. 



Complicating most analyses is the fact that the in- 
ferred dynamical mass is degenerate with the anisotropy of 
the velocity dispersion tensor and the latter is poorly con- 
strained by available line-of-sight velocity data. While this 
degeneracy leaves the full mass profile underconstrained in 
a sta ndard Jeans analysis l|Strigari et al . 2006; Walkc r~et al.l 
2009), its relative weakness near the halflight radius of the 
stellar tracer makes estimates M(R e ) oc nR e a 2 (where R e 
and a are luminous effective radius and global velocity dis- 
persion, respectively, and k is a constant) robust to various 
forms o f anisotropy and/or even to the sha pe of the mass 
profile (|Walker et al.ll2009l ; IWolf et al.ll201fj| ). 

The presence of at least two chemo-dy namically dis- 
tinct stellar subpopulations in several dSphs l|Tolstov et al.l 
1 20041 ; IBattaglia et al.l l201lf ) then provides a unique oppor- 
tunity to measure the slopes of dSph mass profiles, F = 
AlogM/Alogr, direct ly by estimating M(R e ) at two dif- 
ferent effective radii. IWalker fc Penarrubial l|201ll . 'WPll' 
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hereafter) introduce a statistical method that uses estimates 
of stellar positions, velocities and metallicities to estimate 
R e and a for each of two stellar subpopulations within the 
Fornax and Sculptor dSphs, obtaining T — 



2.611^ and 



2.95Z 



respectively. Taken at face value, these mea- 



surements exc lude, with significance ~ 96% an d ~ 99%, re- 
spectively, the lNavarro. Frenk fc White! l| 19971 . 'NFW here- 
after; r ^ 2 at all radii) profile that is often invoked to 
characterise density profiles of cold dark matter (CDM) ha- 
los formed in dissipationless cosmological simulations. WP11 
tested their method against spherical dynamical models with 
various degrees of anisotropy and found that mass estima- 
tors of the form M(R e ) oc R e a 2 systematically overesti- 
mate the enclosed mass more strongly for tracers that are 
more deeply embedded (i.e., more concentrated) in their 
host haloes. This bias implies that slopes T = logM/Alogr 
tend to be systematically underestimated, such that WPll's 
claim of their quoted levels of NFW exlcusion were conser- 
vative. 

However, despite the assumption of spherical sym- 
metry that is common to most dSph studies (excep- 
tions include the axisym metric Schwarzschild analyses of 
Ijardel fc Gebhardj 2012 and the flattened models consid- 
ered bv lAgnello fc Evansll2012h . the composite stellar popu- 
lations of real dSphs are clearly not spherical. The Milky 
Way's 'classical' d Sph satellites have ellipticitie s ranging 
from 0.1 < e < 0.6 (|lrwin fc Hatzidimitrioulll995h . Further- 
more, haloes fo rmed in CDM cosmological simulations ten d 
to be triaxial |Allgood et all l200fj ; IVera-Ciro et al.ll201l[) . 
Therefore, insofar as CDM represents the null hypothesis 
regarding cosmological structure formation, the relevance 
of inferences drawn from spherically-symmetric analyses de- 
pends critically on their robustness to axisymmetric and tri- 
axial cases. 

Here we test the slope measurements of WP11 for ro- 
bustness against non-spherical symmetry. We exploit the 
fact that in a triaxial potential, the energy is an integral 
of the motion and thus we can construct isotropic stellar 
distribution functions (DF) of the form f(E) even within 
triaxial N-body dark matter haloes. W e use the prescrip- 
tion presented bv lLaporte et all i|2013h to build stellar DFs 
with various degrees of concentration within cosmologica l 
CDM haloes produced in the A quarius (|Springel et al. ! 2008) 
and Phoenix (|Gao et al.ll2012T l runs to cover a wide range of 
triaxiality parameters from Milky Way to cluster type envi- 
ronments. Section 2 discusses the numerical simulations and 
method used to generate DFs. Section 3 describes our use 
of samples from these DFs (projected along various lines of 
sight) to examine systematic errors inherent to the WP11 
method for various forms of the tracer number density pro- 
files. We discuss results and conclude in section 4. 



2 NUMERICAL METHODS 
2.1 Dark matter haloes 

For the modelling of dSp h dark matter haloes , we use the 
Aquarius simulations fsee lSpringel et all |2008l ) for details). 
This is a set of six high-resolution dark matter only sim- 
ulations of the formation of Milky Way mass dark matter 
haloes in ACDM. In the level-2 resolution the particle mass 



is ~ 1O 4 M0 and the softening length is e = 65 pc comoving. 
We extract a number of dark matter haloes in the mass range 
10 9 - lO ^fo^Mo, where h = 0.73, using the subhalo finder 
SUBFIND (|Springel et al.l feOQl'). The shape of the Aquarius 
subhaloes have axis ratios which increase with radius and 
which are mildly triaxial with axis ratios < b/a >~ 0.75 and 
< c/a >~ 0.6 at 1 kpc (Vera-Ciro private communication). 
We also complemen t our sample wit h subhaloes drawn from 
cluster simulations (|Gao et al.ll2012T ) to bracket the range of 
possible triaxiality parameters for subhaloes in CDM and 
rescaled masses by a factor of 1000. 



2.2 Generating Tracers 

The weightin g sch eme used here was developed by 



Laporte et al. (|2013l ) and is a generalisation of that of 



Bullock fc Johnston! J2005l ) to triaxial systems. For details 
' Laporte et all (|2013l 7 In short, in order to generate a lu- 



minous stellar profile, we take each simulation particles of 
energy E — hv 2 + $ to simultaneously represent dark mat- 
ter and stars in diferent amounts through the weight func- 
tion (jj(E) — ~ ^E)g(E) ' wnere N is the differential 
energy distribution, g is the density of states and asterisks 
denote stellar quantities. One generates f*(E) through spec- 
ifying the target number density profile v — u(r) and using 
the Eddington formula with an additional approximation 
to deal with the multivalued behaviour of $ = &(r) with 
spherical radius. In this way, the method creates a stellar 
profile which retains contours of the flattening of the total 
potential. 

Figure 1 displays projected number densities and line 
of sight velocity dispersion profiles obtained by sampling 
random projections of DFs (for two different stellar num- 
ber density profiles and different concentrations) calculated 
using the machinery described above. We check for equilib- 
rium by tracking the stellar DFs forward in time for a period 
of 150 Myr, during which the stellar density profiles remain 
stable. 

WP11 originally tested their method using models in 
which stellar populatio ns trace da r k mat ter pot entials char - 
acterised by generalised lHernquistl (|l990l . see also lZha3l996l ) 
profiles: 



v*{r) = v 



and 



Pdm(r) = po 



( 7 »-/3»)/c«» 



(1) 



V r dm J L \ r dm 



•(2) 



For the stellar number densities, WP11 considered Plum- 
mer profiles, («»,/?*, 7*) = (2, 5, 0), which provide good fits 
to dSph surface brightness profile s (jlrwin fc Hatzidimitrioul 
ll995l : lMcConnachie fc Irw in 2006) and has the virtue of de- 
pending on a single parameter, the projected halflight radius 
R e . They also considered alternative profiles that retain a 
luminous core (7* = 0) but fall off more slowly /quickly at 
large radius than do Plummer profile^], with (a*,/?*, 7*) = 



They actually considered models with 7* = 0.1 because mod- 
els 7 = have f(E) < in some regions. In the simulations 
the resolution limit already prevents this from happening for our 
models. 
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(«.,/S.,7.) = (4,4.5,l) - 




1.0 
R [kpc/h] 



Figure 1. Projected number density and line of sigh ve- 
locity dispersion profiles for one Aquarius subhalo using the 
present weighting scheme. Top left: Plummer profiles with r* = 
0.3,0.5, 1.0/i _1 kpc, Top right: (a*,/3»,7») = (4,4.5,1) profiles 
with r* = 0.3,0.5, 1.0ft - x kpc. The vertical dotted line marks the 
point where r = 2.8e. 



(2, 4, 0) and (2, 6, 0), respectively. IStrigari et al.1 l|2010t l have 
shown that cuspy tracer number density profiles provide a 
good match to the observed surface brightness and veloc- 
ity dispersion profiles of the composite stellar populations 
in dSphs. At a fixed half-light radius, a cuspy tracer compo- 
nent would have a lower velocity dispersion than would its 
cored counterpart. In order to test for sensitivity to the in- 
ner profile of the tracer components, here we consider models 
with stellar cusps (a*, /?», 7*) = (4,4.5,1) as well as cored 
Plummer profiles with (2, 5, 0). 



3.1 The bias in the WP mass-estimator: 
systematics 

After calculating DFs as described above, we project each 
model along 100 random lines of sight uniformly sampled on 
a sphere. For each projection angle, we then calculate the 
half-light radius R e of each population. In order to mimic 
the WP11 method, we estimate R E by x 2 -fitting a Plummer 
profile to the tracers. The mass enclosed within R e is then 



M(R e ) oc R, 



V D 2 

— — oc R e a , 



(3) 



where Wi are the N-body particle weights, and a is the global 
velocity dispersion of the tracers. The slope is then calcu- 
latedasr=^iMi^£2). 

In order to check whether the WP11 method continues 
to give conservative limits, Figure 2 displays distributions of 
the bias E\T] = r est — r true over all randomly-chosen view- 
ing angles. In nearly all cases the estimated slope is smaller 
than the true slope, such that the estimated slopes continue 
to represent conservative lower limits. This behavior holds 
regardless of the degree of triaxiality and/or whether the 
light profile is cusped or cored. 

We emphasise that the velocity dispersion that enters 
the WP11 mass estimator, M oc R e a 2 refers to the global 
dispers ion of all stars in a give n stellar sub-population. Re- 
cently, iKowalczvk et all (|2012l . 'K12' hereafter) have found 
that use of a different mass estimator - one that refers to the 
velocity dispersion only of stars inside the half-light radius - 
would give less reliable limits on V, particularly when triax- 
iality is present. We confirm this result using our own DFs 
(Figure 3): indeed, when velocity dispersions are estimated 
using only stars inside R e of their respective subpopulation, 
the estimated slopes have large scatter about the true values 
and do not constitute reliable lower limits. 

In addition to stellar number densities used in studying 
dwarf spheroidals, we also show results from additional tests 
for which we adopted the Jaffe profile, (a*,/?*, 7*) = (2,0,2), 
which can be used to model ellipticals and which has steep 
stellar cusp (7, = 2). In this case, we determine the half- 
light radius by fitting a de Vaucouleurs profile. For these 
cases with steep stellar cusps, we find that that WPll's 
method becomes unreliable when the stellar populations are 
highly concentrated (top panel in Fig. 2); however, for suffi- 
ciently extended stellar populations the method still recovers 
a conservative (i.e., biased towards low values) estimate of 
the slope of the underlying mass profile, albeit with a more 
prominent tail towards positive values. 



MASS MODELLING: 
METHOD 



MULTI-COMPONENT 



The presence of multiple stellar populations in some dSphs 
enables the observer to estimate enclosed masses at two dif- 
ferent half-light radii in the same potential. Testing their 
method on DFs drawn from spherically-symmetric models 
with cored light profiles, WP11 find that masses tend to be 
over-estimated more strongly for more-concentrated stellar 
populations. As a result, the slope T tends to be underesti- 
mated, providing conservative lower limits on the true slope. 
We now use our models /* (E) to test whether this behavior 
holds for the case of triaxial haloes and/or when the tracer 
number density profiles are cusped instead of cored. 



3.2 Why triaxality does not matter so much? 

We can understand the relative insensitivity of the WP11 
method to triaxiality by considering the coupling of esti- 
mated quantities R e and a with respect to projection angle. 
Let us rotate an individual halo in the frame of its body 
axes (as evaluated at a radius of lft _1 kpc) such that the 
major axis lies on the a;-axis and the minor axis lies on the 
y-axis. We then observe it along different polar angles in the 
x — y plane and estimate the half-light radius and velocity 
dispersion via the same x 2 fitting procedure used above. We 
notice that when the velocity dispersion is large (along the 
major axis) the estimated value of R e is at its minimum 
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Figure 2. The stacked bias distribution in the slope determi- 
nation of live N-body dark matter haloes (all observed through 
100 different random lines of sights). The vertical line marks the 
point where F ES t — Vtrue is zero. Left hand: Results for Aquar- 
ius subhaloes for Plummer and (4,4.5,1) profiles (in black and 
dotted magenta respectively) Right Hand:: Phoenix rescaled sub- 
haloes for Plummer, (4, 4.5, 1) and Jaffe profiles (in black, dotted 
magenta and dashed blue respectively). The half-light radii of 
the stellar populations are determined through fitting a Plum- 
mer profile to the number density profile (as assumed in WP11). 
(rel,r e2 ,r eS ) = {0.3,0.5,l-0)h' 1 kpc. 

value and vice versa (Figure 4). This anti-correlation of R e 
and a tends to cancel the effects of triaxiality on the mass 
estimator. Therefore the slope F will be less sensitive be- 
cause at a fixed angle 6 any bias in M(R e ) will cancel out 
in the estimate of T. This is why the biases we recover in 
Figure 2 are similar to those found by WP11 for spherically 
symmetric models. The fluctuations in the mass estimates 
due to triaxiality vary from 10 to 20 percent depending on 
the embededness of the tracer population. 



4 DISCUSSION AND CONCLUSIONS 

We have presented families of isotropic distribution func- 
tions of the form f(E) in triaxial potentials extracted from 
dark-matter-only simulations. These span a range of dark 
m atter density profiles fo r which we have tested the method 
of IWalker fc Penarrubial |201lf ). Our tests show that the 
method is generally able to place conservative limits on 



Figure 3. Same as figure 2 but showing results derived using 
the method used in Kowalczyk ct al. 2012. Clearly losing the 
kinematic information outside the half-light radius of the tracer 
makes the estimator highly unreliable. This method is not that 
used in WP11. 

slopes of mass profiles, even when the light profil es have 
NFW-like cusps as advocated bv lStrigari et al.l l)2010l ). Thus, 
we conclude that triaxiality has little impact on published 
analyses of dSph stellar kinematics that assume spherical 
symmetry. The reason is that R e and a are anti-correlated 
over the range of projection angles, effectively cancelling the 
effects of triaxiality. However, we have found that the WP11 
method can break down if the stellar tracers are highly con- 
centrated and have steeply cusped number density profiles, 
e.g., the Jaffe profiles examined in Section 3.1. Some of 
the haloes which were identified by SUBFIND are strongly 
stripped and the tracer may not be entirely in equilibrium. 
However, using those models, we were still able to recover 
successful limits on the slope of the dark matter density 
profiles. This suggests that tidal stripping does not unduly 
impact the results of WP11. 

Recently, a similar study on the same subject has been 
carried out by K12. Our work differs in three aspects: 

(i) We consider haloes which form within a ACDM cosmo- 
logical context. K12 have considered spherical models which 
get tidally stirred under a static potential. 

(ii) Our models do not have rotation. Many galaxies in 
K12 still retain rotation, which is not observed in dSphs. 

(iii) K12 do not test the robustness of WP11 to triaxal- 
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Figure 4. Amplitude of the variations in R e (black), <r (red) and 
M(r e ) (blue) for different stellar populations (different line styles) 
for an example halo from Aquarius as a function of angle 8, the 
polar angle in the plane of the major and minor axes of the halo 
evaluated at l/i _1 kpc. The anti-correlation in the behaviour of R e 
and a creates almost a cancellation and weak variations in M(R e ). 
We also see that variations are greater the more embedded the 
stellar system similarly to the mass bias observed in WP11. This 
in turn explains why the slopes estimates are still reliable under 
our f(E) models in triaxial potentials. 



ity, but show that a mass estimator based on the velocity 
dispersion within the half-light radius of a tracer can misin- 
terpret the true value of the slope of the total mass profile. 
We confirm their result in Figure 3. 

Finally, we note that given a density profile, there ex- 
ists many possible velocity dispersion profiles which may be 
consistent with the observed data (allowing for anisotropy). 
However, WP11 showed this is not an issue for their method 
under anisotropic Ossipkov-Merritt models but also those 
with constant anisotropy. Combined with the results of our 
current study, the WP11 method seems to be robust to both 
anisotropy and halo triaxiality. 
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